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RANK DEFICIENCY OF KALMAN ERROR COVARIANCE 
MATRICES IN LINEAR TIME-VARYING SYSTEM WITH 
DETERMINISTIC EVOLUTION 


KARTHIK S. GURUMOORTHYtt, COLIN GRUDZIEN§ , AMIT APTEtt, ALBERTO 
CARRASSU" AND CHRISTOPHER K. R. T. JONES§ 

Abstract. We prove that for linear, discrete, time-varying, deterministic system (perfect model) 
with noisy outputs, the Riccati transformation in the Kalman filter asymptotically bounds the rank 
of the forecast and the analysis error covariance matrices to be less than or equal to the number 
of non-negative Lyapunov exponents of the system. Further, the support of these error covariance 
matrices is shown to be confined to the space spanned by the unstable-neutral backward Lyapunov 
vectors, providing the theoretical justification for the methodology of the algorithms that perform 
assimilation only in the unstable-neutral subspace. The equivalent property of the autonomous 
system is investigated as a special case. 
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1. Introduction. The problem of estimating the state of an evolving system 
from an incomplete set of noisy observations is the central theme of the state es¬ 
timation and optimal control theory [7], also referred to as data assimilation (DA) 
in geosciences [6, 20]. In the filtering procedure, based on the concept of recursive 
processing, measurements are utilized sequentially, as they become available [7]. For 
linear dynamics, and when a linear relation exists between measurements and the 
state variables, and when the errors associated to all sources of information are Gaus¬ 
sian, the solution can be expressed via the Kalman filter (KF) equations [8]. The KF 
provides a closed set of equations for the first two moments of the posterior probabil¬ 
ity density function of the system state, conditioned on the observations. In the case 
of nonlinear dynamics, the first order extension of the KF is known as the extended 
Kalman filter (EKF) [7], whereas a Monte Carlo approximation is the basis of a set of 
methods known as Ensemble Kalman filter (EnKF) both of which have been studied 
extensively in geophysical contexts [13, 5]. 

Atmosphere and ocean are example of dissipative chaotic systems. This implies 
the sensitivity to initial condition [11] and the fact that the estimation error strongly 
projects on the unstable manifold of the dynamics [18] which has inspired the devel¬ 
opment of a class of algorithms known as assimilation in the unstable subspace (AUS) 
[23]. In AUS, the span of the leading Lyapunov vectors (to be defined precisely in 
later sections), or a suitable approximation of this span, is used explicitly in the anal¬ 
ysis step: the analysis update is confined to the unstable subspace [16]. AUS has 
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been formalized in the framework of the EKF, EKF-AUS [22], and in the variational 
(smoothing) procedure, 4DVar-AUS [21]. Applications with atmospheric, oceanic, and 
traffic models [24, 3, 17] showed that even in high-dimensional systems, an efficient 
error control is achieved by monitoring only a limited number of unstable directions, 
making AUS a computationally efficient alternative to standard procedures. The AUS 
methodology is based on and at the same time hints at a fundamental observation; 
the span of the estimation error covariance matrices asymptotically (in time) tends 
to the subspace spanned by the unstable-neutral Lyapunov vectors. 

The search for a formal proof of this aforesaid property is the basic motivation of 
the present work which is focused on linear non-autonomous, and linear autonomous 
perfect-model dynamical systems. The main results of the paper are as follows. In 
Theorem 3.5 we show that the error covariance matrices, independent of the initial 
condition, asymptotically become rank deficient in time and then in Theorem 3.7 we 
characterize their null spaces by proving that the restriction of the these matrices onto 
the stable backward Lyapunov vectors converges to zero in time. When restricted to 
the linear, autonomous system with the time invariant propagator A, we establish that 
the stable space of the time independent backward Lyapunov vectors equal the stable 
space of —span of generalized eigen-vectors of corresponding to eigen-values 
less than one in absolute magnitude—in Theorem B.3. Consequently, in Corollary 4.2 
we show that the null space of the error covariance matrices contain the stable space 
of A^ asymptotically. 

The paper is organized as follows. After describing the general notation in Sec¬ 
tion 2, the non-autonomous case is considered in Section 3. The assumptions used 
in proving our main result, other useful results such as the Oseledet theorem, and 
the concepts of observability and controllability for noiseless systems are described in 
Sections 3.1, 3.2 and 3.3. The Theorem 3.5 discussing the rank deficiency of error 
covariance matrices is presented in Section 3.4 and the proof of the Theorem 3.7 using 
the geometric viewpoint of Kalman filtering [2, 25, 1], is detailed in Section 3.5. Sec¬ 
tion 3.6 presents some numerical results buttressing the theorem. Section 4 includes 
the proof of Corollary 4.2 along with a numerical illustration supporting the analytical 
findings for autonomous systems. We conclude in Section 5. 

Although the extension of these results to the general nonlinear case is the object 
of active research [19], the current findings already provide a formal justification to the 
AUS foundation and further motivates its use as a DA strategy in nonlinear chaotic 
dynamics. 

2. Notations. The dimension of the state space is represented by d. For any 
square matrix Z S £,d-y.d gg^ {Ai(^), • • • , Xd{Z)} represent the eigen-values of 

Z where \Xi{Z)\ > ■ ■ ■ > |Ad(^)|. Similarly, let the set {ai{Z), ■ ■ ■ ,ad{Z)} stand for 
the singular-values of Z with <Ji{Z) > ■ ■ ■ > ad{Z). We define the column vectors 
of the matrix Uz = [vi(Z'),- • • ,Vd{Z)] to be the generalized eigen-vectors of Z of 
satisfying the relation ZVz = VzJ{Z) where J{Z) is the Jordan-canonical form of 
Z. In the event that Z is diagonalizable {J{Z) is diagonal), let the entries of the 
diagonal matrix Az = J{Z) symbolize the eigen-values of Z and the columns of Vz — 
the eigen-vectors—be of unit magnitude. Z* denotes the adjoint of Z for the scalar 
product under consideration in and Z^ represents the conjugate transpose of Z. 
For the canonical scalar product (u, v) = uW in Z* = Z^ and when confined to 
the real space where (u,v) = u^v, Z* = Z'^. Unless explicitly stated we assume 
a real vector space endowed with a canonical scalar product. The matrix norm ||Z|| 
we consider is the largest singular-value <Ji{Z) of Z. The notation Z > 0 {Z > 0) 
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is used when Z is symmetric, positive-definite (positive-semidefinite). For any two 
symmetric matrices Z 21 the notation Z\ > Z^ means Z\ — Z 2 > 0. The following 
definitions are useful. 

Definition 2.1 (Real span). The real span of a complex vector w = u -I- iv 
where u, v S is the vector space %v C defined as 

7w = {au -F /3v : a, /3 e K}. 


Definition 2.2 (a-eigenspace). Given an a > 0, the a-eigenspace of a square 
matrix Z denoted hy £°‘{Z) is the real span of the generalized eigenvectors of Z cor¬ 
responding to eigen-values X with |A| < a. 

3. Non-autonomous systems. 

3.1. Set up and Assumptions. We define the general linear non-autonomous 
dynamical system at time n > 0 by 

^n+l — ^n+l^n “f (^-t) 

Yn+l — .^n+l^n+1 T (Jn+l 

where x„ G q„ G p„ G The x„ are the state variables, p„ represents 
model noise, y„ represents observational variables and q„ is the observational noise 
term. The basic random variables {xq, qi, q2, •'' > Pii P2, •'' } are all assumed to be 
independent and Gaussian with 

xo ~ A/'(xo|o, Ao) qn ~ A/'(0, Q„) p„~A/'(0,/), 

such that Aq G is the initial error covariance matrix of the state variable xq, 

Qn G is the observation error covariance matrix at time n, and Fn G 

The matrices Aq, (3„, A„, are known for all time n. Further, A„ and Qn are 

considered to be non-singular, ||A„|| < ca, |lQn|| < cq, and ||iln|| < CH,^n > 1 where 
c^, CQ and ch are positive constants. The model noise error covariance is given by 
Pn = FnFn ■ Unless explicitly stated Aq > 0, i.e, its eigen-values are strictly positive. 

Filtering theory deals with the properties of the conditional distribution, called 
the analysis in the context of DA, of the state x„ at time n conditioned on observations 
Yo-.n = [yij y 2 i ■■•j y™] up to time n where the first observation yi is assumed to occur 
at time n = 1. This conditional distribution provides an optimal state estimate in 
the least squares sense [7]. Under the assumptions of linearity and Gaussianity stated 
above, this conditional distribution is Gaussian, with mean and covariance denoted 
by x„|„ and A„ respectively: 

x„|„ = E[x„ I Torn] , and A„ = E[(x„ - x„|„)(x„ - x„|„)^ | Yom] • 

We also note that the conditional distribution, called the forecast in DA literature, 
of the state x„_|_i conditioned on observations lb:n up to time n is Gaussian with its 
mean and covariance denoted by x„_|_i|„ and respectively: 

Xn+l|n ~ E[Xyj^_i I To:n] , and = E[(Xyi^i X^_|_]^|,.j)(x^^i ^n+l\n'} I ho;n] • 

In this work we concern ourselves with systems that have no model error, i.e, 
Fn = 0 Vn > 1 and investigate the dynamics 

Xn+l — , 


and y^^^i — F{n-\-i^n+i F Qn+i • 


(3.2) 
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We will be interested in asymptotic properties of the conditional error covariances 
and A„. The Kalman filter provides a closed form, iterative formula for obtaining 
these quantities [7]. Under the assumption of no model noise, the update equation 


for the forecast error covariance is 

S„ = AnAn-lA^ . (3-3) 

By dehning the Kalman gain matrix Kn as 

K„ = E„i7j + Q„] , (3.4) 

the analysis error covariance equals 

A„ = (/-iA„i7„)E„. (3.5) 

The update equations for the means are given by 

^n+l|n -^n+l^n\n (3-6) 

^n+l|n+l ^n+l\n 4“ ^n+l {Vn+l Hn+l^n+l\n) ■ (3-7) 

Defining the sequence of matrices Mn as 

Ml = (/ - KiHi)Ai, M„ = (/ - KMA^Mn-i (3.8) 

and writing the propagator i3rn:m+n from time m to time m + n by 

^m-.m+n — ^m+n^m+n—1 ' ’ ’ ^m+li (3-9) 

the analysis covariance at time n can be expressed as 


A„ =(/ - K^H^)A^ ■■■{!- K,H,)A^AoA^ ■ ■ ■ = M„AoB^„. (3.10) 

This equation clearly shows that the asymptotic properties of A„ are closely related 
to those of Bg-n and M„. The notation in equation (3.10) is suggestive of the line 
of argument we will take in the following sections. To outline, we may consider the 
singular value decomposition of the propagator Bq.^ = VnSnU^, and decompose the 
error covariances into a basis of the left singular vectors. In particular, we know that 
this decomposition may be written as a function of the singular values, provided we 
have an appropriate bound on M„ in equation (3.10). Moreover, the left singular vec¬ 
tors of the propagator Bq.^ will become arbitrarily close to the backwards Lyapunov 
vectors of the system. 

The properties of Bq.^ are basically determined by the dynamical system and 
are discussed in the next section, while those of M„ are commonly discussed in the 
context of control theory and are discussed in Section 3.3 where we prove a useful 
bound on its matrix norm in Lemma 3.3. 

3.2. Oseledet’s theorem. Note that the boundedness condition on An implies 
the bound ||i?o:n|| < {cAY^'^n. Then Oseledet’s multiplicative ergodic theorem in [15] 
states that for each non-zero vector u S K'’* the limit 
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exists and assumes up to d distinct values /ii > • • • > which are called the Lyapunov 
exponents. We will assume 


0 > fido+i (3.11) 

so that exactly do < d oi the Lyapunov exponents are non-negative. Further, defining 
the matrices 

^ni^) = , El{m) = [{Bm-.m+n)*Bm-.m+nY" , (3.12) 

Oseledet’s theorem guarantees that the following limits exit, namely 

E^{m) = lim E^{m), (3.13) 

n—>-oo 

E^Tn) = lim £',^(to). (3.14) 

n—>oo 

The eigen-vectors of E^{m) and E^ {'m) represented as the column vectors of L^{m) = 
[li(m),-- - ,1 ))(to)] and [m) = [1{(to),-- - ,l^(m)] respectively are defined as the 
backward and the forward Lyapunov vectors at time m [10]. We note that the 
asymptotic results in later sections will essentially use the backward Lyapunov vectors 
V’ijn). 

The convergence of the individual matrix entries in equations (3.13) and (3.14) 
guarantee the convergence of their characteristic polynomials—whose coefficients are 
well-defined functions of the matrix entries—the roots of which are the eigen-values. 
Therefore, 


lim A^fa lim A^// , 

n-s-oo \ I n-J-oo \ ) 

where we recall that A^ is a diagonal matrix comprised of eigen-values of Z. Using 
the notations from Section 2 we additionally find 

11 Aj {E^irn)) Vj (E^(m)) - U^(m)vj (E^(m)) || < |Aj (EYm)) - Xj {E^{m))\ 

+ \\Eiim)-E\m)\\ 

from which we can infer that 

lim \\Xj {E\m)) Vj (U^(m)) - E\m)vj (U^(m)) || = 0 

leading to lim VEbirn) = = Lf{m). Similarly, lim V^f, . = VEftm) = 

L-f (m). 

Oseledets theorem also asserts the eigen-values of E^{m) or iff (m) do not depend 
on the initial time m, are the same for the forward and backward matrices and relate 
to the Lyapunov exponents as 

fij =\og{Xj{E)), jG{l,---,d} ( 3 . 15 ) 

where we deliberately drop the index m and the superscript 6 or / on E. However, 
the forward and backward Lyapunov vectors are different from each other and they 
also depend on the time to, i.e., L^{k) Y Lf{m) Y L'^(to) Y for fc 7 ^ to. 

Consider the singular-value decomposition Bq.^ = UnSniVn)’^ so that under the 
canonical inner product 

Uf( 0 ) = [{Bo-.nfBo-.n]^ = [VniS^YiVnf]^ = K (^„) " (Un)^ , 
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implying = Vn and 

lim ||vj_„ - 1^(0)|| = 0 (3.16) 

n—>-oo 

where (and similarly Uj^„ below) is the column vector of Vn (respectively t/„). 
Likewise, we obtain 

E'^in) = [Bo..n{Bo..nf]^ = [Un{Snf{Unf]^= Un{Sn)^ {Unf 

from which we can deduce that Vgi, („) = C/„ and 

lim ||uj „ - l5(n)|| = 0. (3.17) 


We also infer that 

{a,{Bo..n))- = XjiE^in)) = Xj{EI{0)). (3.18) 

3.3. Controllability and observability for linear dynamics. The notions of 
observability and controllability are dual notions within filtering problems. Roughly 
observability is the condition that given sufficiently many observations, the initial 
state of the system can be reconstructed by using a finite number of observations. 
Similarly, controllability can be described as the ability to move the system from any 
initial state to a desired state over a finite time interval. Formally stated: 

Definition 3.1. The system (3.1) is defined to be completely observable if 
Vn > 1, 




^ ^ {Bn:n+m) ^n+mQn+m^Ji+m^n:n+m 1 0 


\m—0 


and it is defined to he completely controllable z/Vn > 0, 


( d ^ \ 

det I ^ ( Bn.i^ni-.n+dEn-\-mEn.i-ni {Bn+m-.n+d} j 7^ 0- 

\m=l / 


(3.19) 


(3.20) 


In addition we describe the system as uniformly completely observable (respectively 
uniformly completely controllable^ if equation (3.19) (respectively (3.20)) is bounded 
from zero uniformly in n. 

We will assume that the system in equations (3.2) is uniformly completely observ¬ 
able, i.e., the inequality (3.19) is uniformly bounded away from zero. Note however 
that this system cannot be controllable since the determinant in the equation (3.20) 
is identically zero for a deterministic, perfect-model system as Fn = 0, Vn.The hy¬ 
pothesis of uniform complete observability assures that the error covariance matrices 
remain bounded over time as seen below. 

Lemma 3.2. Suppose that the linear, non-autonomous system (3.2) where the 
initial state Xq has a Gaussian law with mean Xo|o and covariance Ag is uniformly 
completely observable (Definition 3.1). Then the error covariance matrices remain 
bounded for all time, i.e, there exist constants ce and ca such that Vn, ||A„|| < ca 
and ||E„|| < c^. 

Proof. The result is proven for autonomous systems in Kumar [9], Chapter 7, 
equations (2.36) and (2.37). Extension to the non-autonomous case is straightforward 
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by rehashing the steps and changing the constants of the autonomous system to their 
time-varying counterparts. □ 

One should note the recent work of Ni et al. [14] has demonstrated a stronger 
result: in continuous, perfect model systems the assumption of uniform complete ob¬ 
servability is sufficient to demonstrate the stability of the Kalman filter. In particular 
this shows that all solutions to the continuous Riccati equation for any choice of initial 
error covariance are bounded and converge to the same solution asymptotically. This 
strongly suggest the same can be shown for the discrete time system, and we will 
return to this point in our discussion of results in Section 5. 

Utilizing only the boundedness of the error covariance matrices, we demonstrate 
that the matrix stays bounded in the following lemma. 

Lemma 3.3. Consider the uniformly completely observable, perfect-model, linear, 
non-autonomous system (3.2) where the initial state xg has a Gaussian law with 
covariance Ag > 0. Then the matrix Mn defined in equation (3.8) is uniformly 
bounded, i.e., there exist a constant cm such that ||M„|| < CM^^n. 

Proof. We first show that the analysis error covariance matrix satisfies the recur¬ 
sive equation 

A„ = (/ - + KMnKl. (3.21) 


Plugging in the Kalman update equations (3.3) and (3.3), the R.H.S of equation (3.21) 
equals A„ — (A„i7j — KnQn)Kn ■ The equation (4.29) in [4] establishes the equality 
Kn = from which the recursion (3.21) follows; further implying that 

A„ > (/ - KMAnAr^-iAlil - KM'^. 


Recursively applying the above inequality gives A„ > M„AgMj. Decomposing Ag = 
Vao^Ao^Ao employing Lemma 3.2 we find 


MnVAoAl^ 


^ < ||A„|| <CA. 


As ||M„|| < 


MnVAoAl 


^'■Ao ^Ao 


the result follows. Note that as Ag > 0 the 


matrix is well-defined. □ 

Ao 

Bearing this bound in mind we shall proceed to discuss the asymptotic properties of 
the error covariance matrices. 


3.4. The asymptotic rank deficiency of the error covariance. We begin 
by introducing a lemma which allows us to formally describe the collapse of the 
eigenvalues of the error covariance matrix. 

Lemma 3.4. For a given e > 0, let Z & be a symmetric matrix such that 

there is a k < d dimensional subspace W C for which 


sup{||Zu|j : ||u|| = l,u G W} < e. 

Then dim(f'^(Z)) > k where the subspace is in accordance with Definition 2.2. 

Proof. Let {vi, • • • , v^} be an orthonormal eigenvector basis for Z corresponding 
to I Xi{Z) |> • • • >1 Xd{Z) I, and let {ui, • • • , Uk} be a basis for W of unit magnitude, 
such that we write 


d 

Ui=^/3zyVj; ; G {1,2,... ,A:} 
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and the matrix of coefficients 


/ /3i.i 

Pi,2 

Pl,d-k + l 0 


0 \ 

^2,1 

^2,2 

P 2 ,d-k + l P2,d-k+2 


0 


Pk-1,2 ■ 


Pk-l,d-l 

0 

\ Pk.l 

Pk,2 


Pk,d-1 

Pk,dj 


is in column echelon form where for every column index j > d — k + 1, the entries 

= * * * = — 0 

and for every row index I < k, ^ jdf ■ = 1 corresponding ||m/|| = 1. Furthermore, 
as Z is symmetric its eigen-vectors form an orthonormal basis and hence ||Zui||^ = 

d-k+l 

For every I < I < k, setting s = k — I + I we find 


d—k-\-s 

> ||Zu,f = 5] l3l,X^,{Z) > AL,+,(Z) = Xl_i^,{Z). 
i=i 

Hence the k smallest eigen-values in absolute magnitude satisfy 

I Xd{Z) |< ••• <1 Xd-k+i{Z) \<e 

and the result follows. □ 

Theorem 3.5. Consider the uniformly completely observable, perfect-model, lin¬ 
ear, non-autonomous system (3.2) where the initial state Xq has a Gaussian law with 
covariance Ag. Then Ve > 0, 3ni > 0 such that if n > ni, I]„ and A„ will each 
have at least d — do eigen-values which are less than e where d — do is the number of 
negative Lyapunov exponents of the system (3.2), i.e., 

dim (f'^(Sn)) > d — do, and dim (f “^(An)) > d — dg (3.22) 


where the subspace £'^ is in accordance with Definition 2.2. 

Proof. As denoted earlier, let /ii > /i 2 > • • • > be the Lyapunov exponents of 
the system (3.2) where dg < d of them are non-negative. The forward stable Lyapunov 
vectors based at time zero is the set which by definitions (3.13) and 

(3.15) satisfy 


lim — 

n—>-oo 77, 


log ( 


So:„l{ 


( 0 ) 


= Mi- 


(3.23) 


Rewriting the analysis error covariance update equation in terms of the transpose 

A„ = MnlX.oBo,n = Ho:„AoM(f 

we get AnM~'^= Bo-.n and in particular 

A„M-^AgilJ(0) = i?o,4(0). 
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Let us therefore define the sequence of vectors 

(3.24) 

By Lemma 3.3 we know that M„ is bounded above, so that the sequence of vectors 
= M“^Ag ^l^(O) must be bounded below. As such, there is a constant Cw such 
that Cw < ||wj_„||,Vn and j G {do + 1, ..., d}. Choose a p > 0 such that for each 
3 G {do + 1,..., d}, p + < 0. Define = ||^^;"|| . Then for a given e > 0, 3ni 

such that for n > ni 

l|A„w,-„|| = ^^|!Bo:nlJ(0)|| < —< e. (3.25) 

The theorem is therefore an immediate consequence of Lemma 3.4. The proof for 
follows along similar lines. □ 

3.5. Null space characterization and assimilation in the unstable sub¬ 
space. The sequence of subspaces defined by the span of {"v^j,n}'j^do+i 
object of study for the remainder of this section. In particular, we wish to establish 
the connection between this sequence of subspaces and assimilation in the unstable 
subspace which utilizes the backwards Lyapunov vectors. 

Definition 3.6. Define to be the d — do x d — do diagonal matrix with 

diagonal entries given by |Aj (A'^(O)) Also, let us define the following d x 

d — do operators 



(3.26) 

K = [Vdo-Sl,™,- - • ,Vd,„] 

(3.27) 

= Ko+iHv ,ldW] 

(3.28) 

Note that equation (3.17) implies that 

l{m \\U: - = 0. 

(3.29) 


Consider the equation (3.10), namely A„ = MnAoVnSnU^, for the analysis error 
covariance A„ at time n in terms of the matrix M„ and the singular-value decompo¬ 
sition of the propagator Bq-^- Noting that = o'j{Bo.n)vj^n and utilizing the 

relation (3.18) we get 

A„C/^ = M„AoV: {Kf ■ (3.30) 

Likewise, recalling that E„ = A„A„_iA^, we can express the restriction of the fore¬ 
cast error covariances as 

(Kf = A„M„_iAoF„^ (C/*)^ . (3.31) 

Making use of the above relations we now prove one of our main result, which 
states that the norm of the restriction of the analysis and forecast error covariances 
onto the backwards stable Lyapunov subspaces must tend to zero. 

Theorem 3.7. Consider the uniformly completely observable, perfect-model, lin¬ 
ear, non-autonomous system (3.2) where the initial state Xq has a Gaussian law with 
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covariance Aq . The restriction of An and I]„ into the span of the backwards stable 
Lyapunov vectors, tends to zero as n —>■ oo. That is 

lim =0, (3.32) 

n—^oo 

lim =0. (3.33) 


Proof By definition log(Aj(£^'^(0))) = p,j, so that the eigen-values Xj{E^{0)) < 1 
correspond to the stable Lyapunov exponents. Recalling that Xdo+i{El{0)) > ■ ■ ■ > 


Xd{El{0)) we find 


A 


eUo) 


= Xdo+i{El{Q)) and 


lim 

n—>-00 


A 


eI{o) 


= Xdo+iiE^O)) < 1 - 


(3.34) 


Consequent to equation (3.34) we can choose a small 0 < p < 1 and sufficiently large 


ni such that when n> ni, A 1 ~ P- 

The restriction of A„ into the span of the columns of Un is given by the equa¬ 
tion (3.30). Note the column vectors of Vn and Un are orthogonal and of unit norm, 
hence ||R^|| = ||tA^|| = 1- We then find for n> ni 


\\AnUf,iUf,f\\< 


A„/ 


EtiO) 


'||M„||||Ao|| < (1 -p)"cm||Ao||. (3.35) 


Consider, 


WAnLl^ {Lf^f II < ||A„||||L((^ (Ll^f - Wn {Wnf II + WAnWn {Uf^f ||, (3.36) 

and Lemma 3.2 states ||A„|| is bounded. Therefore, 

hm |1 A„l:(^(L^«)^|| =0 (3.37) 


by equations (3.17) and (3.35). This may be similarly stated for the forecast error 
covariance. □ 

The forecast and analysis error covariance matrices for a generic non-autonomous 
system in general do not converge, but the above results entail that asymptotically the 
only relevant directions for the error covariance matrices are the backwards unstable- 
neutral Lyapunov directions validating the central hypothesis made by Trevisan et al. 
[22] in their proposed reduced rank Kalman filtering algorithms. 

An intriguing consequence from equation (3.25) in Theorem 3.5 is the following 
corollary. 

Corollary 3.8. Suppose that for some eg > 0, Ng > 0, and for every 0 < e < eg, 
n > Ng, 


dim (£1'^(A„)) = d — doj (3.38) 

ie: asymptotically the rank deficiency of the analysis error covariance An is exactly 
of dimension d — dg. Then the transformation M“^Ag ^ asymptotically maps the 
forwards stable vectors into the span of the backwards stable vectors 

{liW}j^=do+i asn-^oo. 
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3.6. Numerical results for a 30-dimensional system. Below we provide 
an illustration for this asymptotic rank deficiency property of the error covariance 
matrices. The state space vector x„ and the observation vector y„ have dimension 
d = 30 and g = 10 respectively. This choice is arbitrary and our simulations with 
different d and q have shown qualitatively equivalent results. 

The time-varying, invertible propagators An G the observation error co- 

variance matrices G and the observation matrices C j^ioxso ^11 

randomly generated for sufficiently large n. We employed the Qdi method [10] to 
numerically compute the Lyapunov vectors and the Lyapunov exponents and it was 
found that the number of non-negative Lyapunov exponents was do = 14. Starting 
from a random positive-definite Ag, the sequence (S„, A„) was generated based on 
the Kalman update equations (3.3)-(3.5). For every n we computed the eigen-values 
of An sorted in descending order. 

Figure 3.1 shows the eigen-values of A„ as a function of n. Barring the dominant 
14 eigen-values, the rest converge to zero serving as a visual testament to Theorem 3.5. 
Furthermore, we also calculated the norm |jA„Uj_„||, j e {1, 2, • • • , d} for all n and plot 
them in Figure 3.2. These norm values are unsorted meaning that the topmost line in 
Figure 3.2 represent the values ||A„Ui_„|| and the bottommost line denote ||A„U£;^ji|| 
for different values of n. For j > dg = 14, ||A„Uj_„|| approaches zero suggesting that 
as n —>■ oo, the row space of A„ (and also S„) coincides the space spanned by the 
unstable-neutral, backward Lyapunov vectors, i.e., the bounds in inequalities (3.22) 
are saturated. 


Analysis error covariance eigen-values profile £ Projection onto backward Lyapunov vectors 



Fig. 3.1. Profile of the eigen-values of An. Fig. 3.2. Norm of the projection 
Counting establishes that the bottom 16 eigen-coefficients ||AnUj,ri|| for varying observation 
values converge to zero. time n. 


4. Autonomous linear dynamical systems. 

4.1. Null space characterization for autonomous systems. The noiseless, 
linear autonomous system can be defined from equation (3.2), with the additional 
assumptions that A„ = A, id„ = H, Qn = Q, are fixed matrices for all n — therefore 
the results about the asymptotic rank deficiency property of the error covariance ma¬ 
trices in Section 3 also apply to autonomous systems. However, a stronger statement 
can be made for time invariant systems because the backwards Lyapunov vectors will 
not vary in time. In fact, the result in this section is even valid for the case when only 
the dynamical system is autonomous (A„ = A) but the observation process is time 
dependent (id„ and Qn depend on n). 

Akin to the non-autonomous case we define 

El = [A^{A^)*]^ , El = [{A’^yA'^]^ (4.1) 

and the similarity with equation (3.12) can readily be seen by setting Bm-.m+n = 
fim in equation (3.9) (hence the omission of the time index m). As before, the 
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existence of the limits 


= lim El, Ef = lim E^ (4.2) 

n—>-oo n—>-oo 

is guaranteed by Oseledets theorem [10]. The eigen-vectors of E^ and E^ are called 
the backward and forward Lyapunov vectors, represented here as the columns vectors 
of and L-^ ordered left to right from the most unstable direction—corresponding to 
the largest Lyapunov exponent—to the most stable direction—corresponding to the 
smallest Lyapunov exponent. Specifically, the Lyapunov vectors are defined globally 
and have no dependence on the time in the linear, autonomous case. Without the 
time dependence on the backwards stable Lyapunov vectors, we obtain a stronger 
statement about the asymptotic null space of the covariance matrices. 

Definition 4.1. Let = [l^g+i, • • • ,1^ ■ Note that Theorem B.3 proved 

in Appendix B states tells us that the span of the columns of is equal to 

Corollary 4.2. Consider the uniformly completely observable, perfect-model, 
linear, autonomous system defined from equation (3.2) where = A, but and 
Qn may depend on n and the initial state Xq has a Gaussian law with covariance Aq . 
Then the restriction of the analysis and forecast error covariances onto tend 

to zero as n —)■ 00 . That is 


lim ||A„L'’^(l'^)^1 | =0 

(4.3) 

n—^oo ^ 

lim =0 

n—¥oo 

(4.4) 


Proof. Combining Theorem 3.7 with Theorem B.3 this is a straightforward con¬ 
sequence. □ 

In all our numerical simulations with arbitrary (and completely observable) choices 
of A, H and Q we have additionally observed convergence of A„ and to a fixed 
A and S respectively and seen their null spaces contain as stated by Corol¬ 

lary 4.2 (refer Section 4.2). Considering the recent work of Ni et al. [14], this strongly 
suggests that the classical result of the stable Riccati equation for completely observ¬ 
able and controllable, discrete autonomous systems [9] has an analogue in the case of 
completely observable, perfect model systems. 

4.2. Numerical results for linear autonomous system. We choose a non¬ 
singular matrix A G ]^30x30 consisting of random entries and set do = 12 

of its eigen-values to be greater or equal to one in absolute magnitude. We ran the 
Kalman filtering system long enough and observed that the analysis error covariance 
do converge to a fixed A and then projected A onto the generalized eigen-space of 
. Figure 4.1 plots the absolute magnitude of eigen-values of A sorted in descending 
order (j Ai(A) |> • • • >| Xd{A) j) in blue color and shows the Lyapunov exponents for 
this system in red shade where we note that the number of non-negative Lyapunov 
exponents is exactly 12 tantamount to the number of eigen-values of A greater than 
or equal to one in magnitude. Additionally, it can be verified that the Lyapunov 
exponents are just the logarithm (to the base e) of the absolute magnitude eigen¬ 
values of A. Recalling the definition of the Lyapunov exponents from equation (3.15), 
this equality also lends credence to our Theorem A.3. The plot in Figure 4.2 displays 
j|A (vj (A”^)) jj;j e {1,2,--- ,d} where Vj {AJ') is the generalized eigen-vector of 
Aj(A). Observe that when j > 12, the norm of the projected coefficients is zero 
rendering a visual confirmation to Corollary 4.2. 
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Fig. 4.1. Lyapunov exponents in blue and 4 3 , j^^rm of the projection coefficients 

the magnitude of the eigen-values of A in red. generalized eigen-space of . 


5. Discussion. We have shown that under sequential Kalman filtering, the error 
covariance for a linear, perfect-model, conditionally Gaussian systems asymptotically 
collapses to the subspaces spanned by the backwards unstable Lyapunov vectors. This 
has been known to practitioners in the forecasting community [1], but had yet to be 
stated in precise mathematical terms. In particular, this foundational work validates 
the underlying assumptions and methodology of AUS. 

At the same time, these results open many new questions for ongoing research 
related to AUS algorithms. For instance, the present results do not formally show 
the equivalence of a fully reduced-rank algorithm such as EKF-AUS applied in such a 
setting. The conditions that imply the convergence of the covariance matrices, given 
arbitrary low rank symmetric matrices chosen as initial conditions have yet to be 
established. Recent work strongly suggests that filter stability for discrete, perfect 
model systems can be demonstrated under sufficient observability hypotheses alone 
[14], Determining the necessary hypotheses for stability of the discrete Kalman filter 
with low rank initializations of the prior covariance matrix in perfect model systems 
will be the subject of the sequel to our work. 

Additionally there are conceptual issues to be resolved in bridging the results for 
linear systems to non-linear settings; the former having the advantage of Lyapunov 
vectors being defined globally in space, whereas the formulation must change in a 
non-linear setting, respecting the dependence on the underlying path. Both of these 
directions of inquiry open rich areas for mathematical research and future algorithm 
design. 

While the ultimate goal of DA is a precise estimate of state for chaotic dynamics, 
it is critical to understand the uncertainty of the prediction. An exact calculation of 
the posterior distribution of states for a high dimensional, complex system is compu¬ 
tationally intractable; as computational resources increase, so will model complexity 
and thus computational efficiency alone will not resolve this issue. This work provides 
an idealized, but general framework for future investigations into low dimensional 
approximations for uncertainty calculation. We hope that a precise mathematical 
framework for understanding the nature of uncertainty for linear systems will lead to 
innovative research to surmount these challenges. 

Appendix A. Eigen-values, singular-values and Lyapunov exponents 
of linear autonomous systems. The results established in this Appendix and 
Appendix B should be treated as an independent body of work elucidating the rela¬ 
tionship between various concepts in linear, autonomous systems and not restricted 
to the domain of DA and filtering theory. While these relationships are known and 
can be retrieved from multiple sources in the literature, we have explicitly proved 
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them here for completeness. Readers familiar with these mathematical connections 
may choose to skip through these sections without any loss of continuity. 

Based on the definition of the matrix in equation (4.1) we find \j{El) = 
As —)■ Ef we also have 

lim Xj{El) = lim [crj(A”)]" = Xj{Ef) j € {1,2,-- • ,d} (A.l) 

n—>-oo n—>-oo 

where the eigen-values Xj and singular-values aj are ordered descending in norm. 
Dropping the label for brevity let J = V^^AVa (instead of J(A)) be the Jordan 
canonical form of A. It is straightforward to see that A" = VaJ^V^^ for any integer 
n. The following inequality stated in Theorem 9 of [12] is quite useful. For any two 
square matrices Zi and Z 2 we have 

aj{Zi)ad{Z2) < CTj-(Z 1 Z 2 ) < aj{Zi)ai{Z2) (A.2) 

Since the singular-values of both the matrix and its transpose are the same, it follows 
that 


ad{Zi)aj{Z2) < aj{ZiZ2) < cri(Z'i)crj (Z 2 ). 


Lemma A.l. For any square matrix Z = VzJ{Z)V 2 


-1 


lim [a,{Z-)]^ = lim h(J(Z)-)]^ 


(A.3) 


Proof. Inequalities (A.2) and (A.3) leads to 

ad ( 14 ) ad (14-') < a^iZ^) < a^ (I4) a,{J{Zr). 


Raising each term to the power 1/n and letting n —> 00 proves the result. □ 

Corollary A.2. For any matrix A let E^ be defined as in equation (4-2) and J 
he the Jordan canonical form of A. Then Xj{E^) = lim K(^”)]",je{i,2,--- 4}- 

n—¥oc) 

Proof. The results follows immediately when we employ Lemma A.l setting Z = A 
in conjunction with equation (A.l). □ 

The theorem below establishes the relation between the eigen-values of the time 
invariant propagator A and the limit matrix E^. 

Theorem A.3. For any matrix A let the matrix E^ he defined as in equa¬ 
tion (4-2). Then the eigen-values of E^ equal the absolute magnitude eigen-values of 
A, i.e, XfiEf) = \Xj{A)\J e {1,2, ...,d}. 

Proof. We consider two different cases. 

case 1: A is diagonalizable. When J is diagonal then aj{J) = |Aj(J)| = |Aj(A)|. 

Recalling that Xj{J'^) = [Aj(J)]",Vn, we get [crj(J”)]" = |Aj(A)| and the result 
follows from Corollary A.2. 

case 2: A is not diagonalizable. Let J\{A) denote the Jordan-block of size k x k 
corresponding to an eigen-value A of A of the form 


/A 1 0 • • • 0\ 
0 A 1 ••• 0 


Ja(A) = 


0 0 0 

Vo 0 0 


A 

0 


1 

V 


(A.4) 
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The following lemma is useful in proving Theorem A.3. 

Lemma A.4. For any matrix A let J\{A) be a Jordan block corresponding to 
eigen-value A of A as defined in equation (A.4-)- Then the singular-values of J\{A) 
respect the following equality, namely 


lim K (Jj^)]-=|A| j e ,fc}, 


(A.5) 


i.e, the limiting singular-values are the absolute magnitude of their respective eigen¬ 
values. 

Proof. Following the standard proof technique for equality results we individually 
show that 


Im [aj (J”)]" < |A| 

n—^oo 

and 

lim [ctj (Jj^)]- > |A| 



.. ,fc} 

(A. 6 ) 


• ,fc}. 

(A.7) 


Let the Nilponent matrix N = J\ — XI with = 0. When n> k — 1 we get 


= (A/ -h Nf = 


k-1 

E 

r=0 




Further, the highest singular-value ai{N^) = 1 for r G {0,1, • • • , fc — 1}. If A = 0 then 
J” = 0 when n > fc — 1 and the result is trivially true. Suppose A 0 define ^ = y • 
Using the identity that for any two matrices Z\ and Zi, ai{Zi-\-Z 2 ) < <Ti{Zi)-\-ai{Z 2 ) 
as stated in Theorem 6 of [12], we have 


w)<iAr 


■fc-i 

E 

r—O 


i<5r 


Let |(5| = ef for any 0 < e < |5|. Then 


(j^) < iArc' 


< lArc'^ 


'k-l 


E 

.r =0 

n 

E 

.r—O 


= iAr^'=(i + er 


Raising to the power 1 /n and taking the limit we get 

lim [cTi (Jj^)]" < |A|(l-he). 


(A. 8 ) 


The above inequality is also true for the rest of the singular-values as o'i(.) is the 
largest. Since e is arbitrary the first inequality (A. 6 ) follows. If A = 0 we get the 
desired, stronger equality result in equation (A.5) as the singular-values by definition 
are non-negative. It suffices to focus on the case A 7 ^ 0 where J\ is invertible. 

To establish the reverse inequality (A.7), let Tx be the Jordan canonical form of 
given by 


/a 1 0 ... o\ 

0 i 1 ... 0 

0 0 0 ^1 
Vo 0 0 0 a; 
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Lemma A.l entails that 


lim 


((A‘)’ 


n—>-oo . 

Applying the inequality (A.6) on T\ gives us 


= lim [<^jiTx)]' 


lim [aj (Tj;*)] " < |xy j G {1, 2, • • • , /c}. 


In particular, 


lim 

n—>-oo . 


CTl 


((A')") 


= lim 


< — 


where the equality stems from the fact that for any invertible matrix Z of size k x k 


h {z-^) = 


^k—j + l i^Z) 


We then get 


lim [crfe(Jr)]" > 1^1- 


(A.9) 


Since ak{-) is the smallest singular-value the inequality (A.9) is also valid for the rest. 

□ 

Now to prove Theorem A.3 note that for any n 


(Jl 


J" = 


0 


0 J? 


V 0 0 


0 

jn J 

"'a,/ 


is a block diagonal matrix and the eigen-(singular) values of J” equals the disjoint 
union of eigen-(singular) values of individual Jordan blocks • In accor¬ 

dance with Corollary A.2 and Lemma A.4 we find Vj G {1,2,..., d}, 


X,iEf) = lim K(J”)]^ = |A,(J)| = |A,(A)|. 


□ 

Appendix B. Eigen-spaces and Lyapunov vectors of linear autonomous 
systems. By a suitable coordinate transformation, namely z„ = V({~^x„, studying the 
dynamics x„+i = Ax„ is tantamount to investigating z„_|_i = Jz„ where J = Vj[^AVA 
is the Jordan canonical form of A. Indeed, 


Zn-t-l — dZn — Va AVaV^ X-n — Va ^n+1- 

Corresponding to the definitions of the matrices and E-^ in equations (4.1)-(4.2), 
let Gn = [(J”)* and let G = lim G„. 

n —>-<50 

We consider the two systems in the different d dimensional spaces and Cj 
where the underlying propagators are A and J respectively. Note that as the matrix 
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Va might be complex (though A is real) the dynamics for the propagator J is examined 
in a complex state space. 

Lemma B.l. If the scalar product in Cj is the canonical one namely, (u,v)j = 
u^v, then Vq = Id where Id is the d x d identity matrix. 

Proof. We find it convenient to handle the following scenarios separately, 
case 1: A is diagonalizable. J is diagonal and so is J". In the canonical inner 
product setting the entries of the diagonal are the absolute magnitude entries of 
J. It follows that G is diagonal and Vq = Vj = Id- 

case 2: A is not diagonalizable. As before, consider the Jordan block J\ given 

in equation (A.4) of size k x k corresponding to the eigen-value A. Define G\ = 
1 

lim [{Jx)*J\] ■ Since Gx is symmetric it is diagonalizable and by Theorem A.3 we 

n—yoo 

have Xj {Gx) = |A|, Vj S { 1 , 2 , • • • , k}. As all the eigen-values of Gx are equal, it is a 
scalar matrix and therefore we can choose = Ik- Since 



fGx, 

0 • 

■ ° \ 


0 

Gx, ■ 

0 

G = 





1 0 

0 • 

• GxJ 


the result follows. □ 

Lemma B.2. Under the definition of the scalar products (u,v)j = u^V^Vav in 
Cj and (u, v)^ = u^v in Vq = V^^VEf. 

Proof. For the aforesaid considerations of the scalar products in Cj and 

J* = (vI^Va'^ J^vI^Va and A* = A^ respectively. Recalling that J = V^^AVa we 
have 


(J”)* 


(vivAy'viiA-)^ {vx 


VX^A^)^ VaVX^A^Va 

1 

2n 



As {EI)'^^ = {A^)^ A^ is symmetric, it is diagonalizable by an orthonormal matrix 
V^s and carries a representation = V^,/(A^/)^"’V^j. We find Kc„ = A^^ and 

Vg„ = V)C^I4 ;^,Vn and the result follows by letting n —>■ oo. □ 

Recall the real span 7w from Definition 2.1 bearing in mind the complex gen¬ 
eralized eigenvectors of any matrix Z always occur in conjugate pairs {w,w} with 
7w = TV- We have the following theorem, namely 

Theorem B.3 (Eigenspace equality). For any matrix A let the matrix be 
defined as in equation (4-2). Then for any a > 0 the corresponding a-eigenspaces of 
E^ and A are the same, i.e, £“ (E^) = £“(A). Equivalently, £°‘ (E^) = £“ (^^)- 
Proof By Theorem A.3 we have Aj(G) = |Aj(J)| = |Aj(A)| = Xj{E^). Recall that 
the eigen-values are ordered with Ai(G) and Xd{G) being the largest and the small¬ 
est respectively. Oseledets theorem states that there exits a sequence of embedded 
subspaces 


OcFdC Fd-i C---CFi=Cj 

such that on the complement of in Fj the growth rate is at most Xj{G) 

[15]. The subspaces Fj can be obtained as the direct sum of the eigenvectors Vj(G) 
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as 


T'j — Vd(G) 0 V(i_i(G) 0 • • • 0 Vj(G) 


where 'VjiG) is the eigenvector of G corresponding to Aj(G). Further, though the 
eigenvectors of G depend on the underlying scalar product in Cj, the embedded 
subspaces Tj and the eigen-values Xj (G) are independent of it [10]. 

Corresponding to the two inner-product definitions in Cj, specifically (u, v)j = 
u^v and (u,v)j = we denote the respective eigenvectors with the super¬ 

script symbols 1 and 2. By Lemma B.l we have Vq = Id = V^^Va and Lemma B.2 
declares that Vq = where Ve/ is computed using the canonical inner product 

in For the given a\et q = argmin^- Xj{G) < a. The invariance of the embedded 
subspace J-q to the underlying scalar product signifies that the real span of the vectors 
{V] 4 vJ(G), • • • , VavJ(G)} equal the real span of the vectors {Vav^(G), • • • , Vav^(G)}. 
As Vj S {1,2 ,-•• ,d}, VAvj(G) = Vj(A) and Vav|(G) = Vj(i5^), the result follows. 
□ 
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